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ABSTRACT 

We lay out and apply methods to use proper motions of individual kinematic tracers for es- 
timating the dynamical mass of star clusters. We first describe a simple projected mass es- 
timator and then develop an approach that evaluates directly the likelihood of the discrete 
kinematic data given the model predictions. Those predictions may come from any dynamical 
modelling approach, and we implement an analytic King model, a spherical isotropic Jeans 
equation model and an axisymmetric, anisotropic Jeans equation model. This maximum like- 
lihood modelling (MLM) provides a framework for a model-data comparison, and a resulting 
mass estimate, which accounts explicitly for the discrete nature of the data for individual stars, 
the varying error bars for proper motions of differing signal-to-noise, and for data incomplete- 
ness. Each of these two methods are evaluated for their practicality and are shown to provide 
an unbiased and robust estimate of the cluster mass. We apply these approaches to the enig- 
matic globular cluster omega Centauri, combining the proper motion from van Leeuwen et al 
(2000) with improved photometric cluster membership probabilities. We show that all mass 
estimates based on spherical isotropic models yield (4.55±0.1) x 10 6 Mq[D/5.5±0.2A;pc] 3 , 
where our modelling allows us to show how the statistical precision of this estimate improves 
as more proper motion data of lower signal-to-noise are included. MLM predictions, based on 
an anisotropic axisymmetric Jeans model, indicate for u> Cen that the inclusion of anisotropies 
is not important for the mass estimates, but that accounting for the flattening is: flattened mod- 
els imply (4.05±0.1) x 10 6 M Q [Z?/5.5±0.2fcpc] 3 , 10% lower than when restricting the anal- 
ysis to a spherical model. The best current distance estimates imply an additional uncertainty 
in the mass estimate of 12%. 

Key words: Stellar dynamics - celestial mechanics - Galaxy: globular clusters: individual: 
NGC5139 



1 INTRODUCTION 

Estimating masses of self-gravitating systems, such as star clus- 
ters, robustly and accurately has been of long-standing importance 
in astronomy. Such estimates must be based on the spatial distri- 
bution of some tracer population (e.g. stars) and on its kinemat- 
ics, radial velocities in most circumstances. However, with HST, 
with adaptive optics and the advent of the next generation tele- 
scopes, proper motion data have and will become more accurate 
and reliable within an observation time of a few years. Proper mo- 
tions at the Milky Way's centr e has already proven a boon (e.g. 
ISchodel. Merritt & Eckarj2009l and references therein). Combined 
with radial velocities, proper motions provide 5 components of 
the 6-dimensional phase space for kinematic tracers. Distance esti- 
mates are usually not accurate enough to determine the relative 6th 
phase-space coordinate within the system. Such a wealth of infor- 
mation should help to estimate accurately masses of objects in the 



* E-mail addresses: ritchiedsouza@yahoo.com (RDS); rix@mpia.de 
(HWR); 



local group, even if accelerations are unavailable. Apart from the 
masses, such 2D or 3D kinematic data constrain quite directly the 
orbit distribution and can uncover more subtle dynamic processes 
in star clusters. 

Proper motion data for cluster dynamics have been first imple- 
mented in the 70's (e.g. Cudworth 1976) . However the analysis was 
limited by the data quality at the time. iLeonard & Merritq dl989l) 
showed that proper motion data were sufficient to predict the dy- 
namic mass of the cluster. However, there was no complete proper 
motion data sample of a star cluster at the time. 

At the present, data from HST have reached suf- 
ficient quality t o study internal cluster kine matic s as- 

trometrically dMcNamara. Harrison & Anderson! l2003l 

lAnderson & van der Marell 201(1) . Also, diffraction limited 
imaging at the VLT, Keck and other 8 m class telescopes is now a 
reality, and with the next generation of telescopes coming online 
(e.g. LIN C-NIRVANA came ra at the LBT with 23 m resolution 
imaging, iHerbst et alj |2004) . proper motions are becoming a 
broadly competitive dynamical tool. They are not restricted to 
areas with long time-baseline data, such as ui Centauri (e.g. 
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Ivan Leeuwen et al . 2000). In addition, for many clusters, it is much 
more practical to get proper motions for very extensive data sets 
(1000s of stars!) as opposed to radial velocities. Proper motion 
based mass estimates scale with the distance to the cluster as D 3 , 
thereby requiring good distance estimates. Hence, it is crucial to 
understand how best to use proper motion data for setting mass 
limits and to learn to deal with the systematic and random errors 
inherent to this technique. The main objective of this paper is to 
review and develop methods for mass estimation based on proper 
motions alone. We then apply these to determine the mass of uj 
Centauri. 

The Southern globular cluster w Centauri (NGC 5139) 
is the most massive a nd one of the most flattene d globular 
clusters in the galaxy ( Geyer, Hopp & Nelles 1983). Its large 
mass and size sets it apart from the bulk of globular clus- 
ters in the Galaxy. In recent years, it has evoked interest be- 
cause of the recent discovery of its mu ltiple populations which 
implies a complex formati on history dHilker & Richtler l2000l : 
Hughes & Wallerstein l2000h - uj Centauri may well constitute a 
remnant nucleus, i.e., the former cluster of a now disrupte d 
galaxy dGnedin et al.l |2002|; iMizutani. Chiba & Sakamoto] l2003h . 
iMerritt. Mevlan & MavoJ {l997) reported a definite rotation in the 
cluster. [Anderson! d2002h reports that there is very little mass segre- 
gation in the cluster, which means that the cluster is not in equipar- 
tition. 

Van de Ven et al. (2006; hereafter vdV06) have published 
an extensive study of uj C entauri using the proper motion data of 
van Leeuwen et al. (2000) in addition to radial velocity data pub- 
lished bv lReiins et al. I d2006h employing the Schwarzschild's orbit 
superposition method. It is important that we compare our results 
with vdV06, and seek to understand the differences that may arise. 

There has been consider able debate i n rece nt years about uj 
Centauri's estimated distance. iLongmord dl992h estimates it at a 
distance of 5.1 kpc, while a dynamical distance of 4.8±0.3 kpc has 
been e stimated by vdV06. Using RR Lyrae stars, Del Principe et al] 
( 2006) estimated a photometric distance of 5.5 ± 0.15 ± 0.15 kp<£^ 
which is in agreement with the res ults (5.3 ± 0.3 kpc) in ferred by 
the eclipsing binary OGLEGC 17 dThompson et alJuOOll) . For this 
paper we adopt Del Principe's values of 5.5 ± 0.2 kpc. 

There has also been a considerable debate on the mass of uj 
Centauri. Previous esti mates of the mass of uj C entauri have ranged 
from 2.4 x 10 6 M P> dMandushev etai]|l99lh to 7.13 x 10 6 M o 
dRicheretal]l99ll) . vdV06 has estimated that the mass of the clus- 
ter is (2.5 ± 0.3) x 1O 6 M at a dynamical distance of 4.8 kpc, as- 
suming that the cluster is neaiiy-axisymmetri c. A similar e s timate 



based on an edge-on axisymmetric model by Merr itt et al 



sets the mass at 2.9 x 10 M . On the other hand. lMiocchi 
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2010) 



using spherical symmetry with an anisotropic velocity distribution 
to model the cluster, gets a mass estimate of (3.1 ± 0.3) x 1O 6 M0. 
Setting accurate limits to the mass of uj Centauri is essential to dif- 
ferentiate between the various theories of its formation process. 

Various approaches to go from kinematic data to mass esti- 
mates have been proposed (Binney & Tremain 2008, BT08), dif- 
fering greatly in their complexity. In Section[2] we review some of 
the available mass estimators and then develop better tools geared 
specifically towards proper motion data. For the present analysis, 
we restrict ourselves to the self-consistent case, where the kine- 
matic tracer density, z/* (r), is immediately related to the mass den- 



sity, p(r), that generates the gravitational potential. This implies 
that the cluster is sufficiently relaxed and that dark matter does not 
appreciably change the dynamics of the cluster. In Section [3] we 
first describe the proper m otion data set available for uj Centauri 
dvan Leeuwen et al]|2000l , vLOO hereafter) and construct a sample 
of stars. In Section[4] we apply the tools developed in this paper to 
the constructed sample to estimate the mass of the cluster. Further 
in Subsection l4.3) we use an anis otropic axisymmetric Jeans mod- 
elling technique ( Cappellari 2008) to account for the flattening of 
the cluster. We discuss the implications of our results and the mod- 
elling techniques we use in Subsection l4.4l Finally we summarize 
our results in Section|5] 



2 MASS ESTIMATORS 

The virial theorem has been long and widely used as a technically 
simple way to estimate the masses of self-gravitational systems 
from the motions of individual tracers or constituents. However, 
a major drawback of estimators that are based on the virial theorem 
is that they tend t o be biased and formally have an i nfinite vari- 
ance. That led e.g. lHeisler. Tremaine & Ba hcall dl985h to develop 
a practical "projected mass estimator" for self-gravitating systems. 

While for many systems (e.g. extra-galactic ones), ra- 
dial velocities are the only observable kinematic component, 
iLeonard & Merrittl dl989h extended the estimator to include proper 
motions, and applied this to estimate the mass of the open cluster 
M35. 

In this Section, we shall outline two estimators that have many 
conceptual and practical improvements. We then see how they stack 
up against t he Jeans equation modelling, which exploits more infor- 
mation (e.g. lWolf et al.|[2oToL for a recent application). 



2.1 Projected Mass Estimator (PME) 

Projected mass estimators try to exploit the directly observable, 
projected properties of a self-gravitating system, while being ef- 
ficient and unbiased in obtaining mass estimates of finite variance. 
In practise, such estimators have been derived by taking moments 
of the Jeans equation for a spherical system. For example, a pro- 
jected mass estimator for proper motion velocities can be written 
as 
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1 iDel Principe et alj <2006l) state the error correctly in distance modulus, 
but did not translate it correctly into distances. 



where vr and vt are the radial and tangential proper motion ve- 
locities, flj^Jhe_projected distance and M is the total mass of the 
cluster dLeonard & Merrittll 19891 Eg. 19). 

Note that such a projected mass estimator uses both compo- 
nents of the velocity dispersion, thereby directly accounting for any 
anisotropy in the system. 

Since any projected mass estimator involves a line-of-sight 
projection and radius averaging, care must be taken to have a proper 
radial sampling of the tracers of the system. In this estimator, par- 
ticles at larger projected distances are given more weight than par- 
ticles at smaller projected distances, because of the inability to re- 
construct the true radial distance for particles at small projected 
distance. 

Such projected mass estimators become biased when the ve- 
locity measurements have significant errors, 8v. The resulting mass 
over-estimate grows quadratically with the velocity error. However, 
we can modify the projected mass estimator to take into account the 
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velocity errors. Since < v 2 bs >=< v 2 n 
sonable unbiased estimator would be 
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where <5-u 2 = 28v 2 R + 8v^. 

Note that such mass estimators do not provide a "formula" for 
estimating the error in the inferred mass. However, bootstrapping 
methods can be used to obtain such error estimates. 

Correcting statistically for the velocity errors as in the above 
estimator is a statistical process, and may fail for a small number of 
particles. Any over- or under-estimation of the proper motion errors 
would again bias any estimate. 

Note also, that the simple projected mass estimator and the 
virial mass estimator assume radial symmetry. If the cluster shows 
deviation from radial symmetry or is flattened, the mass estimates 
will be biased. We can quantify how much the projected mass esti- 
mator is biased with increasing flatteni ng by constructing N-Body 
realizations of lowered Evans model (Kui iken & Dubinskilll994l) 
with decreasing potential flattening (q<&). We then apply the pro- 
jected mass estimator (Eq. [TJ to synthetic error-free data from a 
face-on view: the bias in the resulting (spherical) mass estimate de- 
pends on the flattening and on the inclination of the cluster as seen 
in Figure[T] Given that q& « 0.9 corresponds to a three times flatter 
mass distribution, the bias is modest (10%). 

It is important to emphasise that the PME is an averag- 
ing process across the entire radial profile of the cluster. This 
poses practical challenges both at the cluster centre, due to image 
crowding and at the outer edge of the c luster due to conta mina- 
tion from physically unre l ated s tars. As iGenzel et alj fcOOOl) and 
Schodel, Merritt & Eckart ( 2009) have pointed out, moment esti- 
mators must be applied cautiously in cases where the data does not 
extend over the entire length of the cluster, or where the observed 
sample is dominated by stars that are intrinsically far from (but pro- 
jected near to) the centre. 

To summarize, the "projected mass estimator" is very easy to 
implement, but has a number of practical and conceptual limita- 
tions. Therefore, we now turn to explore more (spatial) informa- 
tion. 



2.2 Maximum Likelihood Models (MLM) 

Here we outline an approach to estimate cluster masses that is much 
more geared towards the discrete data that is actually available: 
individual proper motion measurements, with errors that are sig- 
nificant w.r.t. the velocity dispersion. Good estimators should use 
high-precision measurements effectively and not overestimate the 
information content of low-precision velocity measurements. This 
can be accomplished in a framework where we directly calculate 
the likelihood of the data given some model predictions, and then 
maximize it w.r.t. different models (e.g. differing in mass). Such an 
approach must be based on families of suitable models that allow to 
predict the components of the proper motion velocity dispersions 
(both radial and tangential) as a function of the projected radius 
(R). This c an be done by either drawing on an analytic model (e.g. 
lKind [l966). or on an orbit based model (vdV06), or on the Jeans 
equation. If assumptions are made about the symmetry (ellipticity) 
and anisotropy of the cluster, the models reduce to a one param- 
eter family, and the data then simply set the mass normalization. 
These model approaches are well established; what is new here is 
the direct data likelihood estimate. 

Let us presume that we can model the projected proper motion 



velocity distribution at each radius as a Maxwellian distribution. 
Thus the probability of a star in the cluster to have a certain proper 
motion velocity v at the projected distance R in the centre-of-mass 
frame assuming no rotation given that the total mass of the cluster 
is M is: 
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where a is the projected velocity dispersion which scales with the 
total mass of the cluster as a(M, R) = y/Wa(R). 

On the other hand, every observation of proper motion veloc- 
ity Vi has a certain error Svi associated with it. The probability of 
actually measuring the velocity Vi at a projected distance R is a 
Gaussian distribution given by: 



Pobs(Vi \ V,Svi) = 



: exp 
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Thus for a cluster of mass M, the likelihood of measuring a 
proper motion velocity Vi, at a projected distance Ri with an ob- 
servational error Svi is, 



Pi(vi\M) = / p obs (vi\v,8vi)p mode i(v\M)Av. 



(5) 



Note that the convolution in Eq. [5] 'down weighs' the impact 
of the data with large errors in discriminating models of different 
mass. 

The log-likelihood of the cluster to have a mass M given N 
observations is: 



£(Vvi \ M) = ^2 \npi(vi \ M), 



(6) 



where the sum is over both radial and tangential proper motion ve- 
locities. 

The corresponding (posterior) probability distribution for the 
model parameter M is maximised for C ma x of the likelihood func- 
tion. Confidence limits {la) are given by the interval for which 
Cmax — C = 0. 5 (c.f. | Lampton. Margon & Bowverll 19761 whose 
statistic S = 2£; lEfstathiou. Ellis & Petersonlll988ir 

Note that the maximum likelihood methods do not directly in- 
volve the density profile of the tracers, v*(r), in the calculation of 
the mass of the cluster. Hence the observable shape and charac- 
teristics of the tracer density profile, v*{r) need not correspond to 
the actual mass density profile of the cluster, p(r). Rather the max- 
imum likelihood estimator only asks what the proper motions of 
individual stars are given that we have this measurement at radius 
Rt and a certain dynamic model. 

The maximum likelihood estimator can be extended 
to axisymmetric models, if the axisymmetric Jeans equa- 
tions are used to model the velocity dispersions and 
the rotation of the cluster {g inney, Davies & Illin gwortr] 

199d: Ivan der Marel. Binnev & Daviesl 

Emsellem. Monnet & BacorJll994l; ICappellarill2008T) . Equat 
can be suitable modified as: 
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where the mean velocity, v and the velocity dispersions depend on 
the p rojected po si tion ( X, Y) and scale with mass appropriately. 

ICappellaril {2008) has proposed a convenient way to 
model the anisotropic axisymmetric Jeans equation by ex- 
tending the Multi-Gaussian Expansion (MGE) method of 
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Figure 1. The role of flattening in mass estimates: the figure shows the ratio between estimated and true mass values using the PME for the N-Body realization 
of a lowered Evans model cluster (Kuijken & Dubinski 1 994) as a function of the flattening parameter of the potential q$ at various levels of inclination. 
The filled squares are the ratio between estimated and true mass values for a face-on view of the cluster (an inclination of 0°). The open triangles and the 
open circles give the ratio between estimated and true mass values when the cluster is inclined at an angle of 50° and 90° (edge-on) respectively. Note that 
flattening in the potential of 0.9 corresponds to stronger flattening in the density q p m 0.7. 



lEmsellem. Monnet & Baco"nl Jl994l) . We derive the necessary for- 
mulas for the proper-motion velocity dispersion and the rotation in 
Appendix A. 



3 lo CENTAURI DATA SETS 

In this Section, we discuss the available data sets on lo Centauri, and 
the pre-selections and pre-processing that is needed to implement 
the concepts of Section[2] 



3.1 Photometric Data Sets 

No two dimensional surface brightness profile of lo Centauri have 
been published so far. The V-band surface brightness data of 
Mevl anl Jl987^ has been used most in the literature. Recently 
Ferraro et aljj200of) have published a mo re accurate surface density 
profile which differs from lMevlanl J 19871) in i ts outer parts. For thi s 
paper, we use the number density profile of iFerraro et al.l J2006I) . 
which ha s been reproduced in Figure [2] lo Centauri is s ignificantly 
flattened jGever, Hopp & Nellesfl983l[white & S hawl 200fj). The 



clus ter has a tidal radi us = 45' dTrager. King & Diorgovskill 19951) . 

iRev et all J2004) have published wide-field and high preci- 
sion BV photometry in the field of lo Centauri for a large num- 
ber of stars for V > 22. Using Ca and Stbmgren by photometry 
( IRev et alJ feoOO). they separate the field stars (for V > 16) using 
the metallicity-sensitive hk index defined as (Ca -b)-{b - y) plotted 
against (b - y). Because even the most metal-rich stars in lo Centauri 
are r elatively metal-d eficient compared with the typical disk field 
stars, lRevetalJ ( l200i) were able to eliminate foreground field stars 
through the use of the hk index. We use this information below to 
assign photometric membership probabilities. 



3.2 Proper Motions 

vLOO have provided a large set of proper motion data for lo Centauri 
drawing on archival photographic image plates. These proper mo- 
tions are based on 100 plates obtained with the Yale-Columbia 66 
cm refractor, covering the epochs 1931-1935 and 1978-1983. Dif- 
ferential proper motions were obtained for 9847 stars for a limiting 
photographic magnitude of 16.0 for the centre of the cluster and 
16.5 for the outer parts. The precisions of the proper motions range 
from ~ 0.1 mas/yr for the brightest to an average of ~ 0.65 mas/yr 
for the faintest stars. These data have already been used previously 
(vLOO, vdV06) to study the internal kinematics of the cluster. 

Other proper motion data sets also exist. Recently, 
iBellini et all d2009l) have made available the first ground-based 
CCD proper motion catalogue of lo Centauri with a 4 year 
baseline using a Wide Field Imager at ESO's 2.2 m telescope. 
For stars within their saturation limit (V > 14.6) to the 
vLOO f aint limit (V ~ 16.5), the estim ated error is ~ 0.75 
mas/yr. lAnderson & van derMarell feo id) have published high- 
quality proper motions using the HST's ACS for 53,382 stars in 
a central R < 2' field limited to stars on or below the sub-giant 
branch and those brighter than m,FA35w = —11, a few magnitudes 
below the turnoff, using a 2.5- to 4- year baseline. Their typical 
proper-motion error is better than 0.1 mas/yr. 

As we are interested in the total mass of lo Centauri (within its 
tidal radius), we fo cus on the large-area, ground-bas ed data, rather 
than HST data (e.g. l Anderson & van der Marell2010l) . For our mass 
estimates of lo Centauri, we restrict ourselves to the proper motion 
data of vLOO. 

As shown in Figure[2] the errors in the proper motion data be- 
comes comparable to the velocity dispersion profile at magnitude 
16, especially for stars in the outer parts of the cluster. Thus rea- 
sonable and robust estimates of the mass can be obtained at around 
magnitude 16. 
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Figur e 2. Observ ational constraints for proper motion based mass estimates of u) Centauri : (a) Left: the normalized surface number density profile of 
Ferraro et al. 1 2006) in logarithmic scale, (b) Centre: radial (squares) and tangential (triangles) proper motion dispersion profiles as a function of projected 
radius of Sample C at a distance of 5.5 kpc (see Section [53J. The solid and dashed lines represent splines fit through the radial and tangential velocity 
dispersions respectively, (c) Right: typical proper-motion errors of Sample C (see Section 3.2) stars used in Section 4 to determine the mass of the cluster as a 
function of apparent magnitude. 



3.3 Sample Selection and Pre-processing of Data for the 
Mass Estimates 

vLOO divide their proper motion observations into 4 classes based 
on the errors in the proper motion determination which depend 
upon crowding and plate quality. We only use stars with well mea- 
sured proper motions (class 0) since their errors appear to be rela- 
tively well behaved (see Fig. 1 1 of vLOO). Class 1 stars are slightly 
disturbed by a neighbour and are localized in the central part of the 
cluster. A comparision of the velocity dispersion curves in Figure 
[3] indicates that class 1 stars have a higher tangential velocity that 
class stars in the central 5' of the cluster. 



3.3.1 Pre-processing of the Data 

Following vLOO, the centre of the cluster was taken to be a = 
201^69065, 5 = -47!47855. To start off, we assume that the dis- 
tance of the cluster is 5.5 kpc. We project the absolute celestial co- 
ordinates a and S into Cartesian coordinates (relative to the cluster 
centre) onto the plane of the sky perpendicu lar to the line-of-sight 
vector through the cluster centre (Konig 1962): 

x — — /odcos 8 sin Aa 

y — /od(sin 8 cos So — cos S sin So cos Aa) (8) 
and 

d = (sin 8 sin So + cos 5 cos So cos Aa) ~ 1 (9) 

where fo is a scaling factor (fo = 1 to have x and y in radians and 
fo = 180/-7T to have x and y in units of degrees), a and So are the 
cluster centre coordinates, and Aa = a — an and AS = S — So. 

We corrected the proper motion data for apparent rotation aris- 
ing purely from projection effects using a canonical distance of 
5.5 kpc and a systematic l.o.s velocity of 232.02 km/s (vdV06). 
Any solid-body rotation cannot be determined from the differen- 
tial proper motions due to the astrometric reduction process em- 
ployed to measure the proper motions in the first place. However 
it can contribute significantly to the kinematics of the cluster. We 
correct for such a residual solid-body rotation with SI = 0.029 
mas/yr/arcmin following vdV06. 



3.3.2 Construction of Samples 

vLOO calculate the membership probability of each star on the basis 
of its velocity and its projected distance. Choosing cluster members 
based merely on these membership probabilities is fraught with er- 
ror, as they are correlated with radial completeness, mean radial 
distance and the limiting magnitude. Figure[4]demonstrates how the 
completeness and velocity dispersion of the sample correlates with 
the assigned membership probability, thereby affecting the mass 
estimate. Samples with higher cut-off membership probabilities are 
not only radially incomplete but also have lower velocity dispersion 
curves and hence lower estimated masses. 

On the other hand, vdV06 select their sample by choosing 
those stars whose proper motion errors were lower than 0.20 mas/yr 
after a proper motion membership cutoff of at least 68 percent, ob- 
taining a limited sample of 2295 stars. 

For the purpose of kinematic studies it is crucial to determine 
cluster membership in a way that is independent of the star's ve- 
locity. In the subsequent analysis, we will therefore base cluster 
membership predominately on photometric properties. 

In the trade-off between purity and size of photometrically se- 
lected samples, we consider two cases. Sample A: we consider only 
the horizontal branch stars in vLOO data, where field contamina- 
tion is apparently negligible. Selecting horizontal branch stars in 
the ranges —0.1 < B — V < 0.4 and 14.4 < V < 16, we get a 
total of 1053 stars. 

Sample B: we correlate the iRev et al.l J2004r) published sam- 
ple of cluster members based on photometry wi th the vLOO data 
(Sample B). The high precision BV photometry o flRevetalJj2004h 
data c ombined with the Ca and Stomgren by photometry dRev et al.l 
2000), efficiently separates cl uster stars from th e typical disk field 
stars. To match the vLOO and IRev et alj 120041) data-sets, we first 
had to align them. We mapped the two data-sets using the IRAF 
CCMAP routine. We then cross-matched the two data-sets in a four 
vector space (RA, DEC, B and B-V). For the magnitude and colour, 
we gave a large leeway of 0.5 and 0.1 respectively. The size of the 
search window of the RA-DEC was 0.0008° (~ 3")- We matched 
a total of 3321 stars. 

However, this photometric approach leads to manifest outliers, 
stars that photometricaly look like tj Centauri members but whose 
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Figure 3. The projected radial (long dashed line) and transverse (dark dot long dashed line) velocity dispersion of stars of Sample C (as defined in Section 
13.31 of class stars only and the projected radial (dotted line) and transverse (dark short dashed line) velocity dispersion of sample of class 1 stars only as a 
function of radial distance. Sample C is defined later in the text. The sample of class 1 stars was chosen in a similar way to Sample C, but using only class 1 
stars. The velocity dispersion of the sample of class 1 stars extends radially upto 6' . Class 1 stars have a much higher tangetial velocity dispersions (dark short 
dashed line) in the central part of the cluster than class stars (dark dot long dashed line). 



velocity is far greater than the presumed escape velocity. Such stars 
will bias any mass estimates. Therefore, it is important that we also 
clean the sample kinematically of such outliers while introducing 
as little bias as possible. To do this, we calculate iteratively the ve- 
locity dispersion and excise stars outside ±3cr, until all the stars are 
within the dynamical bounds. For Sample A, this clipping removes 
< 10% of the stars leading to 968 members. For Sample B, 2.4% 
stars are removed by the 3<r clipping, bringing the numbers of stars 
in the sample to 3241. 

Using the radial velocity data of iReiins et al.l d2006l) . we can 
derive an independent rough estimate of the number of field stars 
in our sample. For the 3<j cleaned Sample B, out of the total of 
1028 matched stars between the two samples, 22 stars were found 
to be field stars. The probability of the number of fields stars in the 
sample is 2.14%. Given the additional data made available by the 
radial velocity data, we can do better. 

Sample C: both Sample A and Sample B suffer from radial 
incompleteness, with Sample B extending partially upto R= 22'. 
We can extend the sample with the help of the r adial velocity data o f 
1589 cluster stars of oj Centauri published bv lReiins et al.l (2006). 
By using only stars of class 0, we have 1276 stars which are 100% 
kinematically-certified cluster members that extend from 3' to 30'. 
Removing the 1028 common stars and the 22 field stars of Sample 
B, we get a total of 3467 stars. After a 3-sigma clipping, we get 
3384 stars (Sample C). Out of these 3335 stars are within 20'. 

Following vdV06, we only use stars within 20' for the ML es- 
timators. The ML estimators are independent of the spatial com- 
pletness of the sample. Possible field stars with high velocities 
> 20' may unnecessarily bias the ML estimator. 

On the other hand, we use the complete Sample C for the PME 
so as not to bias the estimator at the outer edges of the cluster. In 
the PME, the estimated mass involves an averaging process over 



the entire length of the cluster and hence completness is impor- 
tant. However, any sample so constructed is still incomplete near 
the center of the cluster due to crowding effects. This issue was 
raised in Subsection l2.1l and will be partially addressed in Subsec- 
tion l33H 

Therefore in the final analysis we will only use Sample C stars 
with the ML estimators for the cluster mass estimates. Sample A 
stars are limited in number and have large errors associated with 
their proper motion. Similarly Sample C is superior to Sample B 
in that it has a slightly higher number of stars and fewer contam- 
inated field-stars. Nevertheless, we will also implement the PME 
with Sample A and Sample C in Subsection |4.1| to get a feel as to 
how the estimator fares. 

3.3.3 Velocity Dispersion 

The projected radial and tangential velocity dispersion profiles for 
Sample C as a function of projected distance R is shown in Figure 
fusing radial bins from vLOO (Table 6). For comparison, we show 
the projected radial and tangential velocity dispersion profile of the 
sample of 2295 stars used by vdV06 in Figure [5] In general, the 
velocity dispersions are very similar, even though the evidence for 
a tangential anisotropy in the outer parts is weaker in our (cleaner) 
Sample C. It must also be noted that the velocity dispersions at 
larger radii (> 19') are highly uncertain, due to the small number 
of stars in these radial bins. 

3.3.4 Radial Profile 

Implicit in the methodology of Section [2] is the assumptions that 
the stars with measured kinematics fairly trace the overall popu- 
lation that yields the density profile, v*(r) or S»(i?). However, 
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Figure 4. The number of stars per radial bin and the radial and tangential velocity dispersions of three samples as a function of projected radial distance. The 
thick line represents a predominately photometrically determined Sample C (3384 stars). The long-dashed line represent a sample of stars whose member- 
ship probability is greater than 98.5%(2.5<r) (3881 stars). The short-dashed line represent a sample of stars whose membership probability is greater than 
99.3%(3ct) (3201 stars). The estimated mass of the cluster also varies accordingly as 4.58 ± 0.09, 4.48 ± 0.08 and 3.62 ± O.O6M [D/5.5 ± 0.2fcpc] 3 
respectively. The mass is derived using a MLM of a single isotropic Jeans model. 



crowding effects and poorer colour information at the centre of the 
cluster leads to a bias in the radial sampling of stars with proper 
motions. This affects our sample. For example, we find that fainter 
sub-samples of Sample C lie at larger radii; proper motion estimates 
for them are precluded by crowding at smaller radii. 

To model the surface number density, we fit a King den- 
sity pro file to the surface density data (see Figure [2} pub- 
lished bv lFerraro et alj J2006l) . The probability distribution, P(R), 
( J P(R)R 2 = 1) for stars in each radial bin is shown in Figure[6] 



We can see that Sample C must be radially incomplete both near 
the centre of the cluster and at larger radii. Such radially incom- 
plete samples will affect the mass estimates from those estimators 
which use an averaging process like the PME. 

For the PME, we can attempt to recover back the original den- 
sity profile by repeatedly sampling the same stars. While this is 
possible at the centre of the cluster because of the large number of 
available stars in the sample in those radial bins, it is difficult at 
the edge of the cluster (> 18') where the number of stars found in 
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Figure 5. Radial and tangential projected velocity dispersion as a function of projected distance R for the proper motion sample used by vdV06. The open 
squares and closed triangles represent radial and tangential dispersions for vdV06, respectively. The thicker lines represent splines fit through the radial (solid) 
and tangential (dashed) velocity dispersions respectively. The radial and tangential proper motion velocity dispersion profiles of Sample C (lighter solid and 
dashed respectively) fit to the data shown in Figure 1, are also included for comparison. The typical size of the dispersion errors are shown at the bottom of the 
plot. 



Sample C is much less than the numbers derived from the surface 
density profile. While repeatedly sampling the same stars may re- 
cover in some way the density profile, it may skew the velocity dis- 
tribution. For example, repeatedly sampling stars that intrinsically 
belong to the outer edges of the cluster but which are projected to- 
wards the centre of the cluster, may alter the projected velocity dis- 
tribution at the centre of the cluster. Similarly the limited number 
of stars at the outer ends of the cluster can bias the higher velocity 
moments. 



4 MASS ESTIMATORS FOR cj CENTAURI 

In this section, we connect the methodology of Section [2] with the 
data sets of Section[3]to estimate the mass of uj Centauri. We start 
by applying the PME (Eq.0, then move on to ML-models (Section 
12.21 ) based on a King model and the Jeans equation. 



Sample A predicts (4.84±0.19) x 1O 6 M [D/5.5 ± 0.2fcpc] 3 for 
the same magnitude. 

Note, that the inclusion of the error terms leads to an estimator 
that is virtually independent of Mu m (i.e., signal-to-noise). The 
high mass estimate reveals that the sample is radially incomplete 
especially at the center of the cluster. 

4.2 Simple Maximum Likelihood Models 

We now proceed to combine the dynamical models directly with 
ML estimates for the uj Centauri data to constrain the mass of the 

cluster. 

We start off with the classic King model dKindl 1966b frame- 
work for this task, then use the spherical isotropic Jeans equation 
to predict the velocity dispersion curve; in the subsequent section 
we generalise this to anisotropic, axisymmetric Jeans models 



4.1 Projected Mass Estimator 

To apply the PME, we first have to adjust our Sample A and C to 
fairly trace the overall population as a function of radius (Figure 
[6^. We do so by Monte-Carlo drawing members from Sample C 
so that their radial distribution matches the profile of Ferraro et al. 
2006 (solid line Figure[6]l. 

The resulting PMEs are depicted in Figure|7] with the 1-<t er- 
ror bars calculated from bootstrapping. For bright magnitude lim- 
its, we have only high precision measurements, but few stars. The 
larger error bars demonstrate the lack of information due to the lim- 
ited amount of tracers of the system. 

At magnitude 16, the PME estimates the mass of the cluster 
for Sample C as (4.90±0.12) x 10 6 M o [D/5.5 ± 0.2fcpc] 3 , while 



4.2.1 King Models 

In Figure [8] a King profile dKindl 19621) i 
profile of lj Centauri ( lFerraro et al.l20061) c; 
radius, r c , and the concentration, c. We found r 



is fit to number density 
characterized by the core 
= 141.68" and 



c = 1.22. Note that these values are different from those originally 
derived from Ferraro et al. (2006) as the values they derived cor- 
respond to a dynamic King model (Sigurdsson & Phinney 1995). 
Given the King profile, we can derive a velo city dispersion profile 
using the technique outlined in iKind dl966h or by integrating the 
Jeans equaton. Even though the surface density profile provides a 
relatively good fit, the model cannot reproduce the projected ve- 
locity dispersion. If we nonetheless use this King model and the 
likelihoods equation(6] we obtain mass estimates for different mag- 
nitude limits (shown in the right panel of Figure[8j. 
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Figure 6. The probability distribution (P(R) ■ R 2 ) for stars in radial bins of 1 arcmin derived from the surface density profile of Ferraro et al. (2006, solid 
line). The dashed line represents the radial distribution of those stars in Sample C. The corresponding number distribution of stars for Sample C is also given. 
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Figure 7. Results for the simple Projected Mass Estimator (Eq. |5J as a function of the limiting magnitude (B). The error bars (la) are calculated through 
bootstrapping. The closed symbols represent Sample C stars, the open symbols represent Sample A horizontal branch stars. As the Projected Mass Estimator 
in Eq. 2 does acccount for the growing proper motion errors (as a function of apparent magnitude), the mass estimate depends only at the < 8% level on the 
chosen flux limit. Note that the radial distribution of the kinematic tracers depends somewhat on their magnitude, as faint sources are underrepresented in the 
more crowded inner regions. Hence some magnitude dependence may arise from the differing radial sampling. 



For example at Mu m = 16, the mass estimate is (4.57 ± 
0.09) x 1O 6 M [D/5.5±0.2kpc] 3 , with the error derived us- 
ing likelihood methods. To compare the various models with each 
other, we can define a measure of the scatter as the range of mass 
estimates at various limiting magnitudes divided by the mass esti- 
mate at Mn m = 16. The scatter for the King model amounts to 
0.096. 



4.2.2 Spherical Isotropic Jeans Model 

The isotropic Jeans equation for a spherical system (neglecting ro- 
tation) is: 



d a d<S> 



(10) 
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Figur e 8. King model fit to ui Cen data: (a) Left: a King profile (R c = 141.676", c = 1.224; King 1962, Eq. 14) is fit to the number density jFerraro et ai] 
2006) (b) Centre: corresponding projected velocity dispersion profile of the King model as compared with the radial and tangential dispersions (dashed line 
and dot-dashed line respectively) of Sample C stars as defined in Section [Ol (c) Right: the resulting mass estimate of the MLM (Section |2. 2\ as a function 
of limiting magnitude. The right panel shows that with the proper MLM estimator, the mass estimate varies systematically by nearly 10%, when including 
fainter and fainter stars. 



where a is the radial velocity dispersion, v is the observable num- 
ber density distribution, and 'J/ is the potential generated by the 
matter distribution at the considered point. 

We use t he thr ee dimensional number density distribution by 
iMerritt et alj { 19971) derived using non-parametric analysis (left 
panel of Figure [9}. We assume that mass follows light, and that 
the density profile also defines the potential of the cluster. The ve- 
locity dispersion is then predicted by Jeans equation. This is shown 
in the central panel of Figure [9] The velocity dispersion curve is 
then used to predict the mass of the cluster using the MLM. 

As in Section l4.2.1l the error-bars with bright Mu, n are large 
due to lack of information in the limited number of tracers. At mag- 
nitude 16, the MLM with the Jeans model estimates the mass of the 
cluster as (4.54 ± 0.09) x 10 6 A'/ Q [D/5.5 ± 0.2fcpc] 3 for Sample 
C stars. The scatter of this model as defined earlier in Section l4.2.1l 
is 0.099. 

The resulting MLM mass estimates of the cluster for the King 
model and the spherical isotropic Jeans model are consistent with 
each other. 

The scatter in the results is indicative of the poor fit of the 
spherical models to the complex nature of u Centauri. The results 
of the MLM are also not consistent with the results of the PME 
and those of vdV06. This indicates that we are poorly modelling 
the cluster, failing to take into account the apparent rotation and the 
flattening of the cluster. Hence we now turn to more sophisticated 
models, namely the anisotropic axisymmetric Jeans model. 



4.3 ML Estimates based on anisotropic axisymmetric Jeans 
Models 

Anisotropic axisymmetric Jeans modelling (e.g. ICaPDellarill2003) 
is one way to model uj Centauri without the restrictive approxi- 
mations of isotropy and spherical symmetry. Ample evidence ex- 
ists that neither appro ximation holds; ui Centauri's velocity distri- 
butio n is anisotropic ( |Andersorill2002l; I Anderson & van der Marell 
l201fj ; vLOO; vdV06). iGever. Hopp & Nellesl jl983h have showed 
that u) Centauri is one of the most flattened globular clusters in the 

Galaxy. 

In the case of ICappellaril ( l2008h . the semi-isotropic axisym- 



metric Jeans formalism is generalised to include anisotropy on the 
basis of two assumptions: a constant mass-to-light ratio M/L and 
a constant velocity anisotropy parameter (3 Z = 1 — i>| jv\. Given a 
detailed description of the surface brightness in the form of a multi- 
Gaussian expansion (MGE), the model can predict the shape of the 
second velocity moments (V r ms) on the basis of two free parame- 
ters f3 z and the inclination i of the cluster. We calculate the relevant 
formulae for the second velocity moments for the two components 
of the proper motion in Appendix [B] For the first velocity moment 
(V), additional assumptions are required. Obtaining the rotation of 
the cluster from this formalism involves the need of separating the 
second moments into random and streaming rotation around the 
symmetry axis by the introduction of an additional parameter (and 
hence assumpt ions) for the first order equation (in terms of k\ see 
Equation 35 o fCapp ellaril2008h . analogous to the Satoh (1980) ap- 
proach. The advantage of using a multi-Gaussian expansion (MGE) 
is that it allows us to set the rotation and the anisotropy (bk) 
individually for each Gaussian component. 

Such modelling is based on the two-dimensional surface den- 
sity distribution of the stars, which has not been published as 
such. However vdV06 have derived a multi-Gaussian expansion 
(MGE) of the o ne-dimensional V-band surface brightness profile of 
iMevlarj l fl987l) using 8 Gaussians. They converted this into a two- 
dimensional luminosity distribution by assignin g a projected flat- 
tening q' k to each Gaussian using the data of lGever. Hopp & Nellesl 
dl983h . Unfortunately, these eight flattened Gaussians do not re- 
sult in the elliptic profile of lGev er. Hopp & Nellesl dl983l) . which 
increases upto R — 7' and then falls again towards larger dis- 
tances (Figure [Toj left panel). We rederiv ed a multi-G a ussian 
expansion to the number density profile of iFerraro et al. (2006) 
once more using 8 Gaussians. Using a Levenberg-Marquardt al- 
gor ithm, we optimised the flatte ning coefficients q' k to reproduce 
the lGeverTr Jopp & Nellesl dl983h data as shown in Figure [10] (left 
panel). 

We follow vdV06 in using an isophotal position angle of 100° , 
and align our coordinates axes (x',y') with the observed isophotal 
major and minor axis of to Cen. 

An important parameter in an axisymmetric model is the in- 
clination of the cluster (i) which relates intrinsic and projected flat- 
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Figure 9. Spherical isotropic Jeans model fit to w Cen data: (a) Left: the reconstructed 3-dimensional spherical number density bv lMerritt etaD Jl997h using 
non-parametric analysis shown by the solid line, (b) Centre: the corresponding velocity dispersion curve predicted by isotropic spherical Jeans equation 
(Equation llOt is plotted along with the radial and tangential dispersion (dashed and dot-dashed line respectively) of Sample C stars as defined in Section l3~3l 
(c) Right: the resulting mass estimate of the MLM (Section 12.21 as a function of limiting magnitude. We will show below that relaxing the assumption of 
spherical symmetry, yields significantly lower mass estimates. 



R [ar 




10 15 
R [arcmin] 




1 5.5 16 
Magnitude 



Figure 10. Anisotropic axisymmetric Jeans model: Left: the flattening q' of the projected stellar density as a function of radial distance as derived by Geyer, 
Hopp & Nelles (1983). The dotted line shows the flattening as a function of radius for the MGE parametrization of vdV06, while the long dashed line shows 
the flattening after adopting the parameters in Table [T] Centre: the average velocity dispersion along the major and minor axis (open squares and triangles 
respectively) calculated for Sample C in concentric annuli. The long-dashed line and the short-dashed line depicts the dispersion along the major and minor 
axis respectively of the best-fit anisotropic axisymmetric model. The open circles denote the total projected velocity dispersion calculated in circular annuli, 
while the solid line is the corresponding best-fit from our anisotropic axisymmetric model. Right: the mass of the cluster as a function of increasing magnitude 
using the anisotropic axisymmetric model. Note that the mass estimates are considerably lower than in the previous Figures, where spherical symmetry was 
assumed. 



tening. This can only be constrained by including the line-of-sight 
velocity data, which vdV06 have done. They arrived at a joint con- 
straint for the inclination and the dynamical distance under the as- 
sumption of axisymmetry: Dtani = 5.6(+1.9/ — 1.0) kpc. Us- 
ing their dynamical models (an axisymmetric implementation of 
Schwarzschild's orbit superposition method), they constrain further 
the inc lination at i = 50 ± 3 ° an d the dynamical d istance. Even 
though iKrainovic et al. I d2005h andlCa pp eU ari (2008) have pointed 
out a degeneracy between the inclination and the anisotropy of the 
cluster, we adopt the inclination value from vdV06. The combina- 
tion of the distance and the inclination we have adopted is consis- 
tent with the above relation. 

To fix Kfe (see Equation 35 of Cappellari 2008) required for 
the parametrisation of the random and ordered streaming motion 
around the cluster, we must compare the model projected differen- 



tial rotation with those derived from the cluster data. We use the po- 
lar aperture grid for the proper motions of vdV06. To derive the dif- 
fe rential rotation of the cluster we use the non-parametric approach 
of lMerritt etail d 19971) on the binned data, wh ere we minimise the 
"penalised log likelihood" (see Equation 4 of lMerritt et alj[l997h 
and used a 'thin-plate smoothing spline' routine found in the R- 
package 'FIELDS'. We forced the velocities at large radii (> 19') 
to go down to zero, so that the main features at lower radii do not 
get swamped by the thin-plate spline routine. The differential rota- 
tion velocity contours are shown in the upper panel of Figure QT] 
Note that any solid-body rotation of a cluster cannot be determined 
from the differential proper motions. 

Given the inclination (i), we can proceed to derive ku for each 
Gaussian component by comparing the rotation curves with those 
predicted by the model. The resulting values of Ku are given in Ta- 
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Table 1. The parameters of the 8 Gaussians from the MGE-fits of the num- 
ber density profile of oj Cen as found by Ferraro et al. (2006). Similar to 
Table 2 of vdV06, the second column gives the central surface density (in 
LmqPC~ 2 ) of each Gaussian component, the third column the dispersion 
(in arcmin) along the major axis and the fourth column the projected flat- 
tening of each component which reproduces the flattening of Geyer, Hopp 
& Nelles (1983). The fifth column indicates the coefficient K used to deter- 
mine the rotation of each Gaussian component. 



k 




a' 


q' 


K 


1 


1290.195 


0.47557 


1.000000 


0.0 


2 


4662.587 


1.931431 


0.9991714 


0.0 


3 


2637.784 


2.513385 


0.7799464 


0.4 


4 


759.8591 


3.536726 


0.724126 


1.1 


5 


976.0853 


5.403728 


0.8556435 


0.6 


6 


195.4156 


8.983056 


0.9392021 


0.0 


7 


38.40327 


13.93625 


0.9555874 


0.0 


8 


8.387379 


20.98209 


1.0000000 


0.0 



ble[T] It must be noted that the peaks of rotation do not lie along the 
isophotal major axis along the y direction. If uj Centauri were ex- 
actly axisymmetric, the kinematic minor axis would coincide with 
the isophotal minor axis. However, FigurefTTIshows th at this is only 
approx imately true. A similar feature was noted by Merri tt et al.l 
dl997l) for uj Centauri while deriving the rotation curves from the 
line-of-sight velocities. 

The anisotropic Jeans model allows us to set the anisotropy of 
each Gaussian component (bk). However we follow the practical 
approach of setting a global anisotropy parameter b (see Section 
14.3.11 ), which we can derive from the proper motion data after fixing 
the rotation of the cluster. Similar to Equation[7] we can construct 
a maximum likelihood estimator. The probability for measuring a 
given individual proper motion fi x and fi y , given the parameters M 
and b are: 

1 

Pmodel 

(fi x \X,Y,M, b) 



^/2ty a 2 (M, X,Y,b) 



exp 



-(ti x -v x (M,X,Y,b)\ 2 
2a 2 (M,X,Y,b) 



with a corresponding equation for To obtain the global 
anisotropy parameter b, we marginalise over the mass M. We then 
maximise the likelihood for our Sample C, varying only the pa- 
rameter b. Following standard practice, we can parameterise this 
anisotropy with the variable /3» = 1 — t)| /v\ = 1 — 1/6, where 
z is the line-of-sight direction. Using Sample C, if the inclination 
of the cluster is fixed at i = 50 ± 3°, we find #, = 0.03 ± 0.04 
indicating that overall the cluster is close to isotropy. 

We also calculate the 2-dimensional velocity dispersion pro- 
file of Sample C along the major and minor axis using a thin-plate 
spline (see Equation 11 of Merritt et al. 1997) and compare them 
with the model velocity dispersion profiles (see Equations IB 61 and 
IB It . The results are shown in Figure[T2] The dispersions were cal- 
culated using the polar apertures given in Table 3 of vdV06 and 
the standard techniques outlined in Appendix A and then fed to the 
thin-spline routine. Note that the elongation of the iso-dispersion 
velocity contours along the minor axis is not consistent with the 
assumptions of axisymmetry. 

With an estimate of the global anisotropy parameter /3 Z , we 
can now obtain the mass of the cluster using the MLM (Equation[6j 
in combination with EquationQT]for eacn component of the proper 
motion using Sample C. For each star in Sample C, we calculate 



the first and the second velocity moments for both the components 
of the proper motion (see Equations IB9IIB8IIB6I and iBll . The final 
mass estimate using our own derived values of flattening in TableQ] 
is shown in Figure[l0] The good fit of the model is reflected in the 
low scatter in the mass estimates: the scatter of this model as de- 
fined in |42Tl is 0.76. The mass of the cluster using the anisotropic 
Jeans model is (4.05 ± 0.10) x 10 6 M Q [73/5.5 ± 0.2fcpc] 3 at the 
limiting magnitude B = 16.0. 



4.3.1 Effects of Anisotropy, Inclination, Rotation and Flattening 
on the Mass Estimate 

Within such a modelling context, it would be instructive to separate 
the effects due to inclination, rotation, anisotropy and flattening of 
the cluster. 

While the value of the inclination of the cluster is consis- 
tent with the distance (5.5 kpc) we have adopted and the relation 
73 tan i = 5.6(+1.9/ — 1.0) kpc, we do have a range of accept- 
able inclination i = 45.5(+8.5/ — 5.9) degrees. However, with 
the assumptions of axisymmetry, the observed projected flattening 
of the cluster sets further constraints on the inclination of the clus- 
ter. From Eq. 14 o f lCappellaril ( [2008l) , the inclination is constrained 
by the flattening: i > cos _1 (q[.) Vfc. The fourth Gaussian of the 
MGE in TableQ]sets the lower limit of the inclination at i = 43.6° . 
At i = 45.5°, the mass inferred from Sample C is (3.90 ± 0.10) x 
10 6 M Q [73/5.5 ± 0.2fcpc] 3 . At the lower end i = 44°, the mass is 
(3.82±0.10) x 10 6 M Q [73/5.5 ± 0.2fcpc] 3 , while at the higher end 
oft = 54°, the mass is (4.15±0.10) x 10° A7 Q [73/5.5 ± 0.2/fcpc] 3 . 
Hence, the inclination of the cluster affects the mass estimate upto 
6%. 

By neglecting the cluster rotation in the model, the 
mass inferred based on Sample C is (4.22 ± 0.07) x 
10 6 M Q [73/5.5 ± 0.2fcpc] 3 ; if we include rotation, we get (4.05 ± 
0.10) x lO 6 Af [73/5.5 ± 0.2fcpc] 3 . Hence accounting for rotation 
affects the mass estimate for the cluster by only 4.2%. 

To isolate the effect due to anisotropy, we calculate the 'best' 
mass of the cluster for j3 z = and j3 z = 0.1. We get the same mass 
estimate in both cases. Hence we conclude that in such a modelling 
environment, the anisotropy of the cluster does not contribute much 
to the mass of the cluster. This is because the observed anisotropy 
is very small. 

Accounting for the flattening of the cluster, however, is im- 
portant. This is already apparent from the fact that the masses 
inferred via the flattened models are 10% lower than our ear- 
lier estimates based on spherical models. To explore this further, 
we define for each Gaussian component a projected eccentricity 
e'j = (1 — qf) ' 5 . We then study how the estimate of the mass 
of the cluster varies as a function of the flattening parameter of the 
model in terms of the relative projected eccentricity e' model / 'e' data , 
where corresponds to no flattening, while 1 would correspond to 
the full flattening derived in Table 1. This results in a difference of 
nearly 8.3% in the mass estimates as shown in Figure[T3] which in- 
dicates the importance of properly determining the fattening of the 
cluster to determine its mass. Including also the effects of rotation 
along with the flattening profile of Table Q] further brings the mass 
to (4.05 ±0.10) x 10 6 M Q [73/5.5 ± 0.2fcpc] 3 which is consistent 
with the results of vdV06. 

At the dynamical distance of 4.8 kpc adopted by vdV06, 
this corresponds to a mass estimate of (2.69 ± 0.07) x 
1O 6 A7 [73/4.8fcpc] 3 
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Figure 11. The top panels show the projected rotation velocity fields of w Centauri obtained using a two-dimensional spline smoothing routine, where the axes 
x and y represent the isophotal major and minor axes, respectively. Velocity contours are given in the centre-of-mass frame and are spaced by 1 km/s, while 
distances are measured in arc minutes. The interpolation errors should be < 0.6 km/s. The lower two figures represent the rotation velocity fields in the x' and 
y' direction for our best-fit anisotropic axisymmetric Jeans model. 



4.4 Discussion of oj Centauri Mass Estimates 

Our analysis presented here has made it clear that the major dis- 
crepancy in the estimated masses of u Centauri arises primarily 
due to the choice of the modelling process. By taking into account 
both the flattening of the cluster and its inclination using axisym- 
metric models, and to a lesser extent by accounting for the rotation, 
we can derive consistently a mass estimate that has only ~ 3% 
error for a given distance. 

However, even in our most complex modelling scenario, 
which included both anisotropy and axisymmetry, we had to draw 
on vdV06 to get a cluster inclination. Yet at the same time, we ob- 
served significant deviations in the 2-d velocity dispersion map. 
Our modelling has also shown that compared to the flattening, 
the inclination, anisotropy and the rotation, have less impact on 
the mass estimates. More accurate constraints must be put on 
the flattening of the cluster through more accurate two dimen- 
sional surface density profiles. Further, ui Centauri with its com- 
plex formation history is an ideal candidate for more sophisticated 
Schwarzschild modelling as attempted by vdV06. A perhaps more 



immediate step would be to incorporate radial velocity data in the 
MLM approach established here. 

Irrespective of the modelling context, the largest uncertainty 
in our mass estimates based on proper motion data, arises from 
the adopted cluster distance, which scales the estimated mass as 
~ [D/5.5±0.2fcpc] 3 . In this paper, w e have assume a distanc e 
of 5.5 kpc derived from RR Lyrae stars toel Principe et al]|2006h . 
with random and systematic distance errors of only 3% each. This 
estimate is cons istent with the estim ate of the eclipsing binary 
OGLEGC 17 ^Thompson etal]l200lh . but substantially different 
from the dynamical distance of 4.8 kpc estimated by vdV06. Reli- 
able dynamical distances depend on a host of factors - but crucially 
on a proper estimation of the velocity dispersions, which in turn is 
dependent on the proper selectio n of the sample. The estim ate of 
vdV06 which is lower than that of lDel Principe et ail ( l2006r) . could 
have come about if there is an over-estimation of the proper motion 
dispersion due to unrecognized interlopers or an under-estimation 
of the line-of-sight velocity dispersion. Intrinsic to these dispersion 
estimates is also a proper correction for perspecti ve rotation (which 
is distance dependent) and solid-body rotation. iBono et all ( 2008) 
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Figure 12. Top: the total projected velocity dispersions along the photometric major and minor axis respectively calculated using a thin-plate smoothing spline 
(see Equation 11 of Merritt et al. 1997). The dispersion was calculated using the polar apertures given in Table 3 of vdV06 and then smoothened with a 
thin-plate spline routine. The centres of the polar apertures are indicated by a dot. Bottom: the velocity dispersions as predicted by the best-fit anisotropic 
axisymmetric Jeans model. 



report a systematic offset between reported absolute photometric 
distances and dynamical distances in the clusters ui Centauri and 
47 Tuc, with dynamical distance estimates being lower. Therefore, 
relying on the precise and complet ely independent distance esti- 
mate from Del Princi peTet al.l feOOot) appears preferable. The 3.8% 
distance errors translates into a « 12% error in the mass estimate. 

As this paper is geared as a proof-of-concept paper of what is 
possible with proper motion data, we have steered away from di- 
rectly using the line-of-sight velocity data. However, tackling the 
issues mentioned above warrants the use of the same. In such cases 
however, more sophisticated modelling schemes than what we have 
considered in this paper would be required to determine the dynam- 
ical distance of the cluster, which we postpone to a later paper. 

5 CONCLUSIONS 

In this paper, we implement two independent methods to determine 
the mass of a cluster from the proper motion data. Both the Pro- 



jected Mass Estimator (PME) and the Maximum Likelihood mod- 
elling (MLM) provide unbiased and robust estimates of the mass of 
the cluster. However the use of the PME is limited due to a number 
of practical difficulties in sampling the cluster. For the MLM, we 
use an analytic King model, a spherical isotropic Jeans equation 
model and an axisymmetric, anisotropic Jeans equation model. 

We apply these approaches to the extensive ground-based 
proper-motion data-se t of u Cen t auri b y vLOO. Using the high pre- 
cision photometry of iRev etaL I fcoo4 we construct a kinemati- 
cally unbiased member samples, that should have less contamina- 
tion than the samples used in previous analyses. 

Among the three models used in the MLM context, the 
anisotropic axisymmetric Jeans equation with some rotation is by 
far the best model accounting for the photometric and kinematic 
data of uj Centauri, though there exists compelling evidence of de- 
viations away from axisymmetry. Neglecting the flattening and in- 
clination of the cluster by using spherical models overestimates the 
mass of the cluster by 10%. This accounts also for the large spread 
in reported masses of oj Centauri. 
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Figure 13. Cluster flattening and mass estimates: using the flattening profile of the 8 Gaussians in Table[T] we study the effect of the variation in the estimated 
mass of the cluster as a function of increasing flattening. We do not consider the rotation of the cluster. We define the projected eccentricity e'j = (1 — q' 2 ) 05 
for each Gaussian component. The left end of the x-axis indicates the absence of flattening, while the right end indicates a flattening which replicates the 
projected flattening profile in Table[T] The figure shows the estimated mass of the cluster using the complete Sample C. This shows the effect of modelling 
flattening on the estimated mass of the cluster. Notice that the results reduce back to the simple isotropic spherical case in the absence of flattening (see Figure 
[9}. If we were to include rotation to the modelling process, the mass of the cluster is estimated to (4.05 ± 0.10) X 1O 6 M0 [D/5.5 ± 0.2fcpc] 3 (not shown 
in the figure). 



The mass of the cluster is (4.05 ± 0.10) x 
10 6 M© [D/5.5 ± 0.2kpc] 3 . The largest uncertainty in the 
absolute mass estimate is the distance to uj Centauri which enters 
as the third power. With a distance error of 3.8% this amounts to a 
mass error of ~ 12%. 

While our results are mostly consistent with those of vdV06, 
the difference arises due to our adoption of the new RR lyrae dis- 
tance estimates and due to our more accurate reconstruction of the 
flattening and the shape of the MGE parametrization leading to an 
increase in the mass estimate by 7%. 

This estimate has a 2.5% statistical error, but 12% uncertainty 
arising from the distance errors. 

The authors thank Glenn van de Ven for helpful comments on 
the manuscript. The authors also thank the referee for her patience 
and for helpful advice. 



REFERENCES 

Anderson J., 2002, in ASP Conf. Ser. 265: Omega Centauri, A 
Unique Window into Astrophysics, Main-Sequence Observa- 
tions with HST. pp 87-+ 

Anderson J., van der Marel R.P., 2009, ApJ, 710, 1032 

Bahcall J., Tremaine S., 1981, ApJ, 244, 805 

Bellini A., Piotto G., Bedin L.R., Anderson J., Platais I., Momany 
Y., Moretti A., Milone A.P., Ortolani S., 2009, A&A, 493, 959 

Binney J. J., Davies R. L., Illingworth G. D., 1990, ApJ, 361, 78 

Bono G., Stetson P. B., Sanna N., Piersimoni A., Freyhammer L. 
M., Bouzid Y., Buonanno R., Calamida A., Caputo R, Corsi C. 
E., Di Cecco A., Dall'Ora M., Ferraro I., Iannicola G., Monelli 



M., Nonino M., Pulone L., Sterken C, Storm J., Tuvikene T., 
Walker A. R., 2008, ApJ, 686, 87 
Cappellari M, 2008, MNRAS, 390, 71 

Del Principe M., Piersimoni A. M., Storm J., Caputo E, Bono, G., 
Stetson P. B., Castellani M., Buonanno R., Calamida A., Corsi 
C. E., Dall'Ora M., Ferraro I., Freyhammer L. M., Iannicola G., 
Monelli M., Nonino M., Pulone L., Ripepi V., ApJ, 2006, 652, 
362 

Efstathiou G., Ellis R.S., Peterson B.A., 1988, MNRAS, 232, 431 
Emsellem E., Monnet G, Bacon R., 1994, A&A, 285, 723 
Evans N. W„ de Zeeuw P. T., 1994, MNRAS, 271, 202 
Ferraro F. R., Sollima A., Rood R. T., Orgilla L., Pancino E., Bel- 

lazzini M., 2006, ApJ, 638, 433 
Genzel R., Pichon C, Eckart A., Gerhard O.E., Ott T., 2000, 317, 

348 

Geyer E. H., Hopp U., Nelles B„ 1983, A&A, 125, 359 

Gnedin O. Y, Zhao H., Pringle J. E., Fall S. M.; Livio M.; Meylan 

., 2002, ApJ, 568, 23 
Heisler J., Tremaine S., Bahcall J. N., 1985, ApJ, 298, 8 
Herbst T. M., Ragazzoni R., Eckart A., Weigelt G, 2004, SPIE, 

1042 

Hilker M., Richtler T., 2000, A&A, 362, 895 
Hughes J., Wallerstein G., 2000, AJ, 119, 1225 
King I. R., 1962, AJ, 67, 471 
King I. R., 1966, AJ, 71, 64 

Konig A., 1962, in Hiltner W.A., ed, Astronomical Techniques, 

University of Chicago Press, p.471 
Krajnovic D., Cappellari M., Emsellem E., McDermid R. M., de 

Zeeuw P. T., 2005, MNRAS, 357, 1113 
Kuijken K, Dubinski J., 1994, MNRAS, 269, 13 
Lampton M., Margon B., Bowyer S., 1976, ApJ, 208, 177 



16 R. D 'Souza and Hans-Walter Rix 



Leonard P. J. T., Merritt D„ 1989, ApJ, 339, 195 

Longmore A. J., 1992, in IAU Colloq. 139, 21 

Mandushev G., Staneva A., & Spasova N., 1991, A&A, 252, 94 

Mayor M., et al., 1997, AJ, 114, 1087 

McNamara B.J., Harrison T.E., Anderson J., 2003, ApJ, 595, 187 

Merritt D., Meylan G., Mayor M., 1997, AJ, 1 14, 1074 

Meylan G, 1987, A&A, 184, 144 

Miocchi P., 2010, A&A, 514, 52 

Mizutani A., Chiba M., Sakamoto T., 589, 89 

PageT., 1952, ApJ, 116, 63 

Prada R, Vitvitska M., Klypin A., Holtzman JA., Schlegel D.J., 

Grebel E.K., Rix H.-W., Brinkmann J., McKay TA., Csabai I., 

2003, ApJ, 598, 260 
Reijns R.A., Seitzer P., Arnold R., Freeman K. C, T. Ingerson, 

van den Bosch R.C.E., van de Ven G, de Zeeuw P. T., 2006, 

A&A, 119, 1824 
Rey S.-C, Lee Y.-W., Joo J.-M., Walker A., Baird S., 2000, AJ, 

119, 1824 

Rey S.-C, Lee Y.-W., Ree C.H., Joo J.-M., Sohn, Y.-J., Walker 

A.R., 2004, AJ, 127, 958 
Richer H.B., Fahlman G. G, Buonanno R., Fusi Pecci F, Searle 

L., Thompson I. B., 1991, ApJ, 381, 147 
Schodel R, Merritt D., Eckart A., 2009, A&A, 502, 91 
SuntzeffN., Kraft R„ 1996, AJ, 111, 1913 

Thompson I. B., Kaluzny J., Pych W., Burley G, Krzeminski W., 
Paczynski B., Persson S. E., Preston, 2001, AJ, 121, 3089 

Trager S.C., King I. R., Djorgovski S., 1995, AJ, 109, 218 

van de Ven G., van den Bosch R.C.E., Verolme E.K, de Zeeuw 
PT, 2006, A&A, 445,513 

van der Marel R. P., Binney J., Davies R. L., 1990, MNRAS, 245, 
582 

van Leeuwen F, Le Poole R. S., Reijns R. A., Freeman K. C, de 

Zeeuw P. T, 2000, A&A, 360, 472 
White R. E., Shawl S. J., 1987, ApJ, 317, 246 
Wolf Joe, Martinez Gregory D., Bullock James S., Kaplinghat 

Manoj, Geha Maria, Munoz Ricardo R., Simon Joshua D., Avedo 

Frank F, 2010, MNRAS, 406, 1220 
Wojtak R., Lokas E. L., Mamon G.A., Gottlober S., Prada F, 

Moles M., 2007, A&A, 446, 437 



APPENDIX A: VELOCITY DISPERSION 

The projected velocity dispersion a± of the cluster is needed for 
the kinematic removal of outliers, and can be calculated in the con- 
text of a maximum-likelihood analysis (vdV06). For a given cluster 
member star with velocity Vi at a given projected distance R, the 
probability that the cluster has a dispersion a(R) is given by 



Pi(vi\a(R)) = / p bs(vi\ v,Svi)pdist(v\ a(R))dv, 



(Al) 



where p a b a is the probability of measuring the velocity Vi as defined 
in Equation [4] and pdist is the probability of a star in the cluster to 
have a certain proper motion velocity v given a velocity distribution 
parameterized by the velocity dispersion a(R). 

The projected velocity dispersion cr(i?) is then derived by 
maximising the likelihood function for a given number of stars 
drawn from a concentric ring centred around a projected distance 
R: 



£(Vui| a(R)) = ^2 ln Pi {vi\ a{R)), 



(A2) 



with the la uncertainty resulting from C = 0.5. Usually, Gaussians 
are a good low-order approximation for the velocity distribution 
Pdist{v), thus simplifying the process. 

We can also modify the distribution function to account for 
th e probability of the star to belon g to the cluster (first introduced 
bv |Pradall2003l;IWoitak et aT]|2007h : 



Pdist(v) = ap v + [1 - a] 



2v n 



(A3) 



where a is the probability of the star to belong to the cluster, p v 
is the regular Gaussian distribution and v max is the maximum pro- 
jected velocity in the particular annular ring where the dispersion 
is being calculated. 



APPENDIX B: ANISOTROPIC JEANS EQUATION 

Following the formalism o f lCappellaril J2008h . we can calculate the 
second moments and the mean velocity along the two per pendicular 
axes projected on the plane of the sky. ICappellaril (2008) have al- 
ready indicated that the projected proper motion disper sion can be 
written via single quadratures. Similar to equation 28 of lCapp ellari 
(2008), we can derive the analogue expressi ons for the proper 
motion s using the formulas in Appendix A of lEvans & de Zeeuwj 
(1994). The y component of the velocity dispersion is: 
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(Bl) 



where we retain the definitions of uok, A, B, C and T>. For the sake 
of completeness we repeat these definitions once again. 
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V = l-b k ql - [(1 - b k ) C + (1 - q 2 ) b k ] u. 

Similarly, the x' component of the velocity dispersion is 
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(B6) 



We can also obtain the proper motion first v elocity moments o f 
the whole MGE model similar to equation 39 of lCappellaril d2008l) . 
The MGE model uses an anisotr opic analogue of S atoh (1980) ap- 
proach, as seen in equation 37 of Cappellari (2008), where k fixes 
the rotation for each Gaussian component. 
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The velocity dispersion can then be calculated as: 
" (BIO) 

(BID 



2 — 2 

a,,/ = u , — iv 



We used the adaptive quadrature routines found in QUAD- 
PACK. Both the integrands in equation [B8| and |B9l can be ill-defined 
for certain combinations of anisotropy and rotation coefficients n k 
towards the centre of the cluster. 



